#ifndef AIRLINE_HPP
#define AIRLINE_HPP
//=========================================================
#include <stdio.h>
#include <math.h>
#include <string.h>
#include <unistd.h>
#include "cube_spline.hpp"
//=========================================================
//  QUANTUM STATE
//=========================================================
struct qstate{
  int num;
  double cx;
  double cy;
  qstate *next;
};
struct speed_len{
  double speed;
  double len;
  speed_len *next;
};
//=========================================================
//  CLASS AIRLINE
//=========================================================
//  ВТ - воздушная трасса
//  КС - квантовое состояние
//=========================================================
class airline{
  char *name; //название ВТ
  char *sidstar; //тип ВТ: SID-STAR
  double length; //длина ВТ
  int quant_state_count; //количество КС
  double quant_state_len; //длина продольного эшелона
  point *first; //набор точек ВТ (список)
  double **array_point; //набор точек (массив) для сплайна
  int point_count;
  qstate *quant_states; //набор центральных точек КС
  speed_len *sl;
public:
  airline(){ //конструктор по умолчанию
    this->name = NULL;
    this->sidstar = NULL;
    this->first = NULL;
    this->quant_states = NULL;
    this->length = 0;
    this->quant_state_count = 0;
    this->quant_state_len = 0;
    this->array_point = NULL;
    this->point_count = 0;
    this->sl = NULL;
  }
  //-------------------------------------------------------
  //  Перевод GPS координат в декартовы
  //  phi-широта lat
  //  psi-долгота lon
  //-------------------------------------------------------
  //градусы.минуты->радианы
  double degmin_to_rad(double angle){
    return (M_PI/180)*angle;
  }
  //-------------------------------------------------------
  //сферическая модель Земли
  point *spher_to_dekart(double phi, double psi){
    point *tmp;
    double R = 3440; //радиус Земли = 6371км или 3440м.миль
    double x = R*cos(phi)*cos(psi);
    double y = R*cos(phi)*sin(psi);
    tmp = new point;
    tmp->next = NULL;
    tmp->x = x;
    tmp->y = y;
    return tmp;
  }
  //-------------------------------------------------------
  //-------------------------------------------------------
  //эллипсоидальная модель Земли
  //WGS84 => CK42
  point *wgs84_sk42(double lat, double lon, double h){
    point *tmp;
    tmp = new point;

    double a_wgs84 = 6378137;
    double b_wgs84 = 6356752.3142;
    double e_2_wgs84 = 6.69437999014*pow(10,-3);
    double f_wgs84 = 1/298.257223563;

    double a_sk42 = 6378245;
    // double b_sk42 = 6356863.019;
    //double e_2_sk42 = (pow(a_sk42,2)-pow(b_sk42,2))/pow(a_sk42,2);
    double f_sk42 = 1/298.3;
    
    //double ro = 206264.8062;
    double delta_x = -25+0.36;
    double delta_y = 141-0.08;
    double delta_z = 80-0.018;
    
    /*double w_x = ;
    double w_y = ;
    double w_z = ;
    double m = 1;

    double delta_a = a_sk42 - a_wgs84;
    double delta_e_2 = e_2_sk42 - e_2_wgs84;
    double a = (a_sk42 + a_wgs84)/2;
    double e_2 = (e_2_sk42 + e_2_wgs84)/2;
*/
    double M = (a_wgs84*(1-e_2_wgs84))/pow((1-e_2_wgs84*pow(sin(lat),2)),3/2);
    double N = a_wgs84/sqrt(1-e_2_wgs84*pow(sin(lat),2));

    //double delta_h = delta_x*cos(lat)*cos(lon)+delta_y*cos(lat)*sin(lon)+delta_z*sin(lat)-a_wgs84*(a_sk42-a_wgs84)/N+(f_sk42-f_wgs84)*(a_wgs84/b_wgs84)*N*pow(sin(lat),2);
    double delta_lat = (1/(M+h))*(-delta_x*sin(lat)*cos(lon) - delta_y*sin(lat)*sin(lon) + delta_z*cos(lat) + ((a_sk42-a_wgs84)/a_wgs84)*N*e_2_wgs84*sin(lat)*cos(lat) + (f_sk42-f_wgs84)*(M*a_wgs84/b_wgs84 + N*b_wgs84/a_wgs84)*sin(lat)*cos(lat));
    double delta_lon = (1/(M+h))*(-delta_x*sin(lon)+delta_y*cos(lon))/cos(lat);

    tmp->x = lat + delta_lat;
    tmp->y = lon + delta_lon;

    return tmp;
  }
  //проекция Гаусса-Крюгера
  point *projection_GK(double lat, double lon, double h){
    point *tmp;
    tmp = new point;
    
    tmp = wgs84_sk42(lat,lon,h);
    lat = tmp->x;
    lon = tmp->y;

    lat = (M_PI/180)*lat;

    double n,l,x,y;
    modf((6+lon)/6,&n);
    l = (lon-(3+6*(n-1)))/57.29577951;
    double sl = sin(lat);
    x = (6367558.4968*lat - sin(2*lat)*(16002.8900 + 66.9607*pow(sl,2) + 0.3515*pow(sl,4) - l*l*(1594561.25 + 5336.535*pow(sl,2) + 26.790*pow(sl,4) + 0.149*pow(sl,6) + l*l*(672483.4 - 811219.9*pow(sl,2) + 5420.0*pow(sl,4) - 10.6*pow(sl,6) + l*l*(278194 - 830174*pow(sl,2) + 572434*pow(sl,4) - 16010*pow(sl,6) + l*l*(109500 - 574700*pow(sl,2) + 863700*pow(sl,4) - 398600*pow(sl,6)))))))/1852;
    y = ((5 + 10*n)*pow(10,5) + l*cos(lat)*(6378245 + 21346.1415*pow(sl,2) + 107.1590*pow(sl,4) + 0.5977*pow(sl,6) + l*l*(1070204.16 - 2136826.66*pow(sl,2) + 17.98*pow(sl,4) - 11.99*pow(sl,6) + l*l*(270806 - 1523417*pow(sl,2) + 1327645*pow(sl,4) - 21701*pow(sl,6) + l*l*(79690 - 866190*pow(sl,2) + 1730360*pow(sl,4) - 945460*pow(sl,6))))))/1852;

    tmp->next = NULL;
    tmp->x = x;
    tmp->y = y;
    return tmp;
  }
  //-------------------------------------------------------
  //добавление координат (х,у) в конец списка
  void add(double lat,double lon,double height,double speed,bool spher = false,bool ellipse = false){
    point *tmp1, *tmp2;
    // printf("speed in air = %g\n", speed);
    if(spher){
      tmp2 = spher_to_dekart(this->degmin_to_rad(lat),this->degmin_to_rad(lon));
      tmp2->speed = speed;
    }else if(ellipse){
      tmp2 = projection_GK(lat,lon,height);
      tmp2->speed = speed;
    }else{
      tmp2 = new point;
      tmp2->next = NULL;
      tmp2->x = lat;
      tmp2->y = lon;
      tmp2->speed = speed;
    }
    tmp1 = this->first;
    if(tmp1 != NULL){
      while(tmp1->next != NULL){
        tmp1 = tmp1->next;
      }
      tmp1->next = tmp2;
    }else{
      this->first = tmp2;
    }
    this->point_count++;
  }
  //-------------------------------------------------------
  //удаление всех точек ВТ
  void del_all_points(){
    point *tmp;
    tmp = this->first;
    while(tmp != NULL){
      tmp = this->first->next;
      delete this->first;
      this->first = tmp;
      this->point_count--;
    }
  }
  //-------------------------------------------------------
  //изменение названия ВТ
  void set_name(const char *value){
    size_t n = strlen(value);
    char *tmp;
    tmp = new char[n];
    for(size_t i=0;i<n;i++){
      tmp[i] = value[i];
    }
    tmp[n] = '\0';
    this->name = tmp;
  }
  //-------------------------------------------------------
  //изменение типа ВТ SID/STAR
  void set_sidstar(const char* value){
    size_t n = strlen(value);
    char *tmp;
    tmp = new char[n];
    for(size_t i=0;i<n;i++){
      tmp[i] = value[i];
    }
    tmp[n] = '\0';
    this->sidstar = tmp;
  }
  //-------------------------------------------------------
  //получить название ВТ
  char *get_name(){ return this->name; }
  //-------------------------------------------------------
  //получить тип ВТ
  char *get_sidstar(){ return this->sidstar; }
  //-------------------------------------------------------
  //получить длину ВТ
  double get_len(){ return this->length; }
  //-------------------------------------------------------
  //получить список КС
  qstate *get_quant_states() {return this->quant_states; }
  //-------------------------------------------------------
  //получить количество КС
  int get_quant_state_count(){ return this->quant_state_count; }
  //-------------------------------------------------------
  //получить длину продольного эшелона
  double get_quant_state_len(){ return this->quant_state_len; }
  //-------------------------------------------------------
  //вычисление длины ВТ
  void calc_len(){
    point *tmp = first;
    double D = 0; //длина ВТ
    if(tmp != NULL){
      while(tmp->next != NULL){
        D += (sqrt(pow(tmp->next->x-tmp->x,2) + pow(tmp->next->y-tmp->y,2)));
        tmp = tmp->next;
      }
    }
    this->length = D;
  }
  //-------------------------------------------------------
  void get_speed_len(double tau = 60){
    speed_len *list_of_speed, *ls;
    point *tmp;
    tmp = this->first;
    double curr_speed = 0;
    list_of_speed = NULL;
    if(tmp != NULL){
      while(tmp != NULL){
        if(curr_speed != tmp->speed){
          ls = new speed_len;
          ls->speed = tmp->speed;
          ls->len = tau*tmp->speed/(60*60);
          ls->next = NULL;
          // printf("ls->speed = %g, ls->len=%g\n",ls->speed,ls->len );
          curr_speed = tmp->speed;
          list_of_speed = this->sl;
          if(list_of_speed != NULL){
            while(list_of_speed->next != NULL){
              list_of_speed = list_of_speed->next;
            }
            list_of_speed->next = ls;
          }else{
            this->sl = ls;
          }
        }
        tmp = tmp->next;
      }
    }
  }
  //-------------------------------------------------------
  void get_speed_quant_state_count(){
    point *cur,*next;
    speed_len *sl = this->sl;
    int count_states = 0;
    int qc = 0;
    double curr_speed = sl->speed;
    double s = 0;
    cur =  this->first;
    printf("quant_state_len = %g\n",sl->len );   
    if(cur != NULL){
      while(cur->next != NULL){
        next = cur->next;
        s += sqrt(pow(next->x - cur->x,2)+pow(next->y - cur->y,2));
        if(s >= (sl->len/2)){
          if(qc%2 == 1){
            this->add_qs(next->x,next->y,count_states);
            count_states++;
          }else{
            if(curr_speed != cur->speed && sl->next != NULL){
              sl = sl->next;
              curr_speed = sl->speed;
              printf("quant_state_len = %g\n",sl->len );   
            }
          }
          qc++;
          s = 0;
        }
        cur = next;
      }
    }
    this->quant_state_count = count_states + 1;

  }
  //-------------------------------------------------------
  //вычисление длины продольного эшелона
  void calc_quant_state_len(double tau = 60){
    //tau секунд время одной операции взлет/посадка, время занятости ВПП
    double speed = 250; //миль в час, скорость ВС (const)
    this->quant_state_len = speed*tau/(60*60);
  }
  //-------------------------------------------------------
  //вычисление количества КС
  void calc_quant_state_count(){
    this->quant_state_count = floor((this->length/this->quant_state_len) + 1);
  }
  //-------------------------------------------------------
  //вычисление количества точек (x,y)
  int count_point(){ return point_count; }
  //-------------------------------------------------------
  //заполнение массива точек (из списка)
  void list_to_array(){
    int n = this->point_count;
    double **array;
    array = new double*[n];
    for(int i=0;i<n;i++){
      array[i] = new double[3];
    }
    point *tmp;
    int i = 0;
    tmp = this->first;
    if(tmp != NULL){
      while(tmp != NULL){
        array[i][0] = tmp->x;
        array[i][1] = tmp->y;
        array[i][2] = tmp->speed;
        // printf("\t точки x=%g  y=%g speed = %g\n",array[i][0],array[i][1],array[i][2]);
        i++;
        tmp = tmp->next;
      }
    }
    //удаление старых данных из массива, если они там были
    if(array_point != NULL){
      int k = round(sizeof(this->array_point)/(3*sizeof(double)));
      for(i=0;i<k;i++){
        delete this->array_point[i];
      }
      delete this->array_point;
    }
    this->array_point = array;
    // for (int i = 0; i < this->point_count; ++i){
    //   printf("-> x=%g\ty=%g\tsp=%g\n",array[i][0],array[i][1],array[i][2] );
    // }
  }
  //-------------------------------------------------------
  //построение кубического сплайна
  void spline(FILE *fi, FILE *flog){
    point *new_points, *tmp;
    list_to_array();
    new_points = cube_spline(this->count_point(), this->array_point, fi, flog);
    this->del_all_points();
    while(new_points != NULL){
      this->add(new_points->x,new_points->y,0,new_points->speed);
      tmp = new_points;
      new_points = new_points->next;
      delete tmp;
    }
  }
  //-------------------------------------------------------
  //получение центральных точек для КС и запись в список состояний
  void set_quant_state_point(){
    point *cur,*next;
    double S = 0;
    int qc = 0;
    int count = 0;
    cur =  this->first;
    if(cur != NULL){
      while(cur->next != NULL){
        next = cur->next;
        S += sqrt(pow(next->x - cur->x,2)+pow(next->y - cur->y,2));
        if(S >= (this->quant_state_len/2)){
          qc++;
          if(qc%2 == 1){
            this->add_qs(next->x,next->y,count);
            count++;
          }
          S = 0;
        }
        cur = next;
      }
    }
  }
  //-------------------------------------------------------
  //добавление координат КС в конец списка
  void add_qs(double x, double y, int n){
    qstate *tmp,*tmp1;
    tmp = new qstate;
    tmp->cx = x;
    tmp->cy = y;
    tmp->num = n;
    tmp->next = NULL;
    tmp1 = this->quant_states;
    if(tmp1 == NULL){
        this->quant_states = tmp;
    }else{
      while(tmp1->next != NULL){
        tmp1 = tmp1->next;
      }
      tmp1->next = tmp;
    }
  }
  //-------------------------------------------------------
  //печать данных о ВТ
  void print(FILE *fd){
    fprintf(fd,"Name: %s\n", this->get_name());
    fprintf(fd,"Type: %s\n", this->get_sidstar());
    point *tmp;
    int count = 0;
    tmp = this->first;
    if(tmp != NULL){
      while(tmp != NULL){
        count++;
        fprintf(fd,"№%i \t %g \t %g \t %g\n",count, tmp->x, tmp->y,tmp->speed);
        tmp = tmp->next;
      }
    }
  }
};
//=========================================================
#endif